Bending rigidity of stiff polyelectrolyte chains: a single chain and a bundle of 

multichains 
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We study the bending rigidity of highly charged stiff polyelectrolytes, for both a single chain 
and many chains forming a bundle. A theory is developed to account for the interplay between 
competitive binding of counterions and charge correlations in softening the polyelectrolyte (PE) 
chains. The presence of even a small concentration of multivalent counterions leads to a dramatic 
reduction in the bending rigidity of the chains that are nominally stiffened by the repulsion between 
their backbone charges. The variation of the bending rigidity as a function of fo, the fraction of 
charged monomers on the chain, does not exhibits simple scaling behavior; it grows with increasing 
fo below a critical value of fo- Beyond the critical value, however, the chain becomes softer as 
fo increases. The bending rigidity also exhibits intriguing dependence on the concentration of 
multivalent counterion U2\ for highly charged PEs, the bending rigidity decreases as ri2 increases 
from zero, while it increases with increasing ri2 beyond a certain value of ri2. When polyelectrolyte 
chains form a iV-loop condensate (e.g., a toroidal bundle formed by N turns (winds) of the chain), 
the inter-loop coupling further softens the condensate, resulting in the bending free energy of the 
condensate that scales as N for large N. 

PACS numbers: 61.20.Qg, 61.25.Hq, 87.15.Da 



INTRODUCTION 



Like-charged polyelectrolytes (PEs) such as DNA and 
F-actin undergo a dramatic compaction to highly or- 
dered structures in the presence of multivalent countc- 
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E Il7l Il8l Ha I20L l2l|. Theoretical description of the 
phases of highly charged PE condensates is complicated 
due to the interplay between the chain deformations and 
counterion-mediated interactions. This complex inter- 
play is often not treated adequately. In the classical work 
on stiff PEs , the effect of the condensed counterions is 
simply to renormalize the backbone charge density. More 
recently, Park et al. [23j and Ubbink and Odijk [24[ have 
considered the bending rigidity of DNA strands that can 
form bundles. The underlying assumption in these cal- 
culations is that total bending energy can be written as a 
sum over all bending energies from individual chains. It 
is, however, not obvious when this is valid, since electro- 
static interactions are long-ranged and the counterion- 
mediated attractions are not pairwise additive. In a pre- 
vious study, we considered charge fluctuation effects on 
the bending rigidity of randomly charged polymers |25j . 
Notably we have shown that charge fluctuations along 
the polymer chains give rise to a reduction in the bend- 
ing rigidity as compared to the cor resp ondin g uncharged 
cases. Several other authors 1261 127L 123. |29| have also 
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considered a similar problem of the bending rigidity of 
a single stiff PE. All the existing approaches, however, 
do not provide much insight into a more realistic case, 
namely bundles of PE chains, since they do not take 
inter-chain coupling into account. Note that the attrac- 
tion that softens PEs also induces interchain collapse. In 
fact, counterion- induced interactions are not pairwise ad- 
ditive and thus the result for a single chain can easily be 
invalidated by the interchain coupling. The key to under- 
standing the bending rigidity of PE condensates lies in 
a consistent treatment of both interchain and intrachain 
couplings. 

In this paper, we present a theoretical model for de- 
scribing the bending rigidity of highly charged stiff PEs, 
both a single chain and many chains which form bundles. 
Throughout this paper, the bending rigidity is estimated 
in terms of the persistence length, which is a length scale 
over which the chain does not change its direction ap- 
preciably. The persistence length is simply the bending 
rigidity divided by fcgT, where ks is the Boltzmann fac- 
tor and T is the temperature. In particular, we study the 
interplay between preferential adsorption of multivalent 
counterions and charge correlations in determining the 
bending rigidity of PE chains. We begin by first consid- 
ering a single chain case. We show that, below a certain 
critical value of the electrostatic strength, the PE chain 
becomes stiffer with the increasing interaction strength. 
Above the critical value, the chain becomes softer as 
the interaction strength increases. Thus a simple scal- 
ing function does not describe the bending rigidity. At 
first glance, this is quite surprising because our results 
imply that the bending rigidity of condensates can be re- 
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duced by increasing the backbone charge density. The 
bending rigidity also varies nonmonotonically with the 
concentration of multivalent counterions, n^. For suf- 
ficiently highly charged PEs, the electrostatic bending 
rigidity drops to a large negative value (signalling the 
collapse of the PE) upon a small increase in 712 from 
zero, while it decreases in magnitude with increasing n% 
beyond a certain value of rii . In other words, these exists 
an optimal value of at which the PE is most efficiently 
softened. In contrast to the previous work l2*fij. we 
consider a bundle of multichains as well. In some cases, a 
long single PE chain can fold back on itself (in a rodlikc 
bundle) or wind around itself (in a toroidal bundle) to 
form a bundle. In the case of rodlikc bundles, there is 
extra complexity arising from hairpin turns of the chain 
at both ends of the bundle, while, in toroidal bundles, 
non-parallel portions of the chain complicate the prob- 
lem [23J- Here we ignore these comlexities; in the re- 
sulting picture, the distinction between rodlike bundles 
and toroidal bundles becomes minor. With this simpli- 
fication, we consider toroidal bundles without a loss of 
generality. Whether a bundle is formed by a single chain 
or many different chains, we can consider the bundle to 
be formed by many loops; each loop in a toroidal bundle 
corresponds to one complete turn (wind) of the chain. 
As a result, all monomers in a loop can be considered as 
interacting with each other through intra-chain (or intra- 
loop) interactions, since any portions of the loop do not 
wind around each other. We then discuss the effect of 
coupling between loops in the bundle on the bending 
rigidity. Our theory implies that, the bending free en- 
ergy of each loop in the bundle cannot be added up to 
give the total bending free energy. We also show that the 
bending free energy cannot be computed by summing up 
the interactions over all pairs of loops. The theory rather 
suggests that the inter-loop coupling further softens the 
condensate. This is due to the many-body nature of elec- 
trostatic interactions induced by counterions. As a result, 
the bending free energy of the condensate varies linearly 
with the number of loops. 



II. MODEL AND BENDING RIGIDITY 
A. general formalism 

Our calculations are based on a model adopted in 
Ref. [H, |2(| . The PEs are uniformly charged chains, 
forming a bundle of ./V loops. Each loop is assumed to 
consist of M monomers of length 6, and each monomer 
carries either a negative charge of — e or 0. In reality, 
different loops can have different numbers of monomers, 
but this polydispersity will not influence the persistence 
length as long as the loop size is somewhat larger than the 
screening length assumed here. If /o is the frac- 

tion of charged monomers, then each monomer is asigned 
a charge — e/o on average. Besides monovalent salts, 
there are (Z : 1) multivalent salts as well in the solution. 



Counterions are divided into two classes, "condensed" 
(i.e., they are bound to the PE chain) and "free" (i.e., 
they move freely in the solution) 0, 0] . The charge on 
loop j of monomer s (in units of electronic charge e) can 
assume the values 



&(s) 



-fl + m\ + m 3 2 Z, 



(1) 



where m^mj = 0,1,2,3,... are, respectively, the num- 
bers of monovalent and multivalent condensed counteri- 
ons at the site. Due to condensed counterions, the charge 
distribution along the chain is heterogeneous. The result- 
ing charge destribution is charac terize d not only by the 
average charge per each site, q~j = qj(s) = — fo+fl+Zfl, 
but also by the variance in the charge per site Aj = 



Note here that 



(. . .) is an average over all realizations of the charge dis- 
tribution, and fl and /| are the average numbers of 
the monovalent and multivalent condensed counterions 
on loop j per site, respectively. The electrostatic inter- 
action Hamiltonian can be expressed in terms of qj(s): 



-«|rj(s)-r./(s') 



„ N M 



(2) 



where Vj (s) describes conformations of PE chains in the 
bundle. The Bjerrum length £b — e 2 / ATreksT is the 
distance at which the electrostatic interaction between 
two elementary charges e is equal to the thermal energy 
fcflT, and e is the dielectric constant of the solution. The 
Debye screening parameter is k 2 = 8it£bI where / is the 
ionic strength. 

The partition function is a sum over all charge vari- 
ables and conformations. We treat the free ions within 
the three dimensional Debye-Hiickel (DH) theory, and 
consider backbone charges together with their condensed 
counterions as forming one-dimensional Debye-Hiickel 
(DH) systems |3, [l0| • Thus our approach should be con- 
trasted with the approach of Nguyen et al. [2|j, which 
is based on a strongly correlated liquid (SCL) model of 
condensed counterions. In Appendix A, we argue that 
the major contribution to the electrostatic bending rigid- 
ity can arise from long-wavelength charge fluctuations as 
captured by our DH approach, even when the condensed 
counterions are strongly correlated. This is because small 
bending is more efficiently felt by long-wavelength fluc- 
tuations, even when counterion ordering is mainly driven 
by short-wavelength fluctuations. However, aspects of 
SCL picture will become more pronounced at low tem- 
peratures. 

For a given conformation of a PE described by [^-(s)], 
we can trace over the charge variables. If we de- 
fine the following block matrix Q[rj(s)] whose element 

Q jfs ArAs)] = A 7 1 %'^+^ ^^ r ^J , (3) 
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then the electrostatic free energy functional becomes 



k B T 



N M 



2 X! J2w'{^j ls jj' s s S > 

jj'=l ss< 

±log{detA,Q[r,( S )]} 



+ 



+self-energy. 



(4) 



Note that the derivation of this result is analogous to 
the derivation of the free energy of bundles of rodlike PE 
chains which have been studied extensively 0, El O, 0, 
l3"H - This result can readily be obtained by performing 
standard Gaussian integrals. So the intermediate steps 
that lead to this result are omitted. The first term in 
eq 0] corresponds to the screened repulsion, while the 
second term comes from charge fluctuations. The last 
term comes from j = j' and s = s', and is a self energy 
that should be subtracted from the free energy in eq0] 
Note here that the self energy does not contribute to the 
persistence length, which measures a free energy cost for 
bending or the free energy of the PE chain with reference 
to a rodlike conformation. This is because the self-energy 
is the same for all chain conformations. 

When PEs are near the rod limit, it is reasonable to 
assume that they are bent uniformly. The central axis of 
PEs can then be parameterized by Tj(s) as $j(s) = s/R, 
where R is the radius of curvature. Note that this 
parametrization amounts to adopting a ground state 
dominance approximation in the chain conformations. 
With this parameterization, the free energy with respect 
to the reference rod conformation can be obtained. To 
this end, we write 

Qjfss' - QjfsAdj{ s )} 

= A- 1 8 jr 5(R$) + 



R 2 $ 2 + R 2 f 



xe 



-k jR 2 tf 2 +R 2 ., 
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-kJr 2 ^ 2 +r 2 ., 

V 33 

R 2 $ 4 kR 2 $ a 
{R 2 d 2 + R 2 ? ,) 3 / 2 + RH 2 + R 2 „, 



,(5) 



where Rjj> is the distance between chains j and /. Note 
that this equation is similar to eq 5.2b in Ref |25| . ex- 
cept that the coupling between loops is included in this 
case. It can be shown that the underlined term in eqEl 
is much smaller than the first two terms near the rod 
limit, i.e., kR ^> 1, by the factor of 1/R 2 . If we de- 
note Q to be the corresponding matrix for the rod- 
like case, i.e., the first two terms in eq [3 then we can 

thus write Q[#(a)] = Q [l+ Qq 1 (QW(s)] - Q ) } 

and consider the term in [•] as an expansion parameter. 
Convergence of the expansion is assured if kR ^> 1. Both 
Q -1 [tf] and detQ[tf] can be computed in powers of 1/R . 



When PE chains form a bundle of N loops, charges on 
one loop correlate not only with others on the same loop 
but also with charges on different loops. In other words, 
the intra-loop and inter-loop charge fluctuations are cou- 
pled with each other. Consistent treatment of intra- and 
inter-loop correlations is hard even when all loops are 
perfectly rigid and parallel (in the absence of conforma- 
tional deformations) [10| . In the following subsection, we 
will provide explicit solutions for the case of N = 1. The 
effects of loop-loop interactions on the bending rigidity 
of PE condensates, namely bundles of PE chains, is dis- 
cussed in Subsection C. The electrostatic contribution to 
the bending free energy of a bundle, i.e., the change in the 
electrostatic free energy due to bending is approximated 
as 



AT, 



clcc 



jj' ss'E-d ^ ^ 

+ll2Y,[(Qo 1 -^)(Q-Qo)] jjss , (6) 

where the subscripts j and j' run over loops. Here, the 
first term comes from the repulsion between charges, 
which is screened by both condensed counterions and 
added salt. The second term, which vanishes as Aj — > 0, 
represents the attraction due to charge fluctuations in 
the monomeric charges. It should be noted here that the 
two terms in eq|S]are put on an equal footing. What is 
suppressed here is non-Gaussian fluctuations which are 
responsible for short-wavelength charge correlations. We 
argue in the appendix that Gaussian fluctuations can eas- 
ily dominate the AJ- e i ec as long as nb is small. This 
justifies our neglect of non-Gausian fluctuations in the 
computation of AT e i ec (thus the persistence length) . 

When kT 1 is somewhat smaller than the loop size, then 
the matrix element {Qq ) j -i ss i m eqEl for example, de- 
pends on \s — s'\ for given j and j' . This enables us 
to further simplify the bending free energy in eq HO by 
summing it with respect to s. To this end, let us define 
matrices D(s) and °D as 



k B T 



1 
12 



dss e 



4 -K Js 2 +R 2 ., 

4 e V jj' 



{s 2 + R 2 3 ,f/ 2 s 2 + R 



Djji(s) = — s 4 e V 33' 
12 



{s 2 + R%,f/ 2 s 2 + R 2 r 



(7) 



and °M and M(k) as 



M jf (k) 



AT 5 



jj' 



AT 1 ^, + 2^K Q {R jf VK 2 + k 2 ), (8) 



where £ = £s/b and Kq(x) is the zeroth-order modi- 
fied Bessel function of the second kind. In terms of 
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°D, D(s), °A4, and Ai(k), the bending free energy can be 
partly diagonalizcd. In other words, the free energy can 
be Fourier transformed from s to its Fourier conjugate k. 
To order 1/R 2 , we obtain 



AK 



lcc 



knT 



L 



(9) 



where 

(hi 



b ^ A,-Ay 

33 



N 

+ E 

x cos ks 



(°M) 



°D (°M) 



dk 



33' 







dsD jf (s) / — [M-W-A], 



2tt 



(10) 



where Ajji = AjSjj/. Despite its apparent complexity, 
the expression in eq 1101 can readily be derived in the 
same spirit as the bundle free energy 0, 0, In 
Appendix B, we outline the detailed intermediate math- 
ematical steps that lead to eq^] Note that the bending 
free energy cannot be written as a sum of two-body in- 
teractions over all pairs of loops. This is consistent with 
the earlier finding that rod-rod interactions in polyelec- 
trolyte solutions are not pairwise additive 0, nil. it 
is not surprising to see that the breakdown of pairwise 
additivity is also manifested in the bending free energy. 



B. single-loop cases: N = 1 

For N = 1, the bending free energy given in ea 1101 can 
be further simplified: 



L 



tQSF 



{CAf 

[l + 2£A^oM] 2 12 
dk 2K {bVk 2 + k 2 ) 
2tt 1 + 2^AK t) {b^/WT^ 2 ) 



ds(s 2 + s 3 k)- 



cos ks. 



(11) 



In this equation and in what follows we drop the subscript 
1 in Ai if it refers to a single-loop case. This expression 
should be compared with the OSF result, i.e., £osf = 
£sq 2 /46 2 k 2 , where eq is the "renormalized" charge per 
site and, in our notation, is q = (—fa + f\ + Zf-i)- In 
OSF theory, the effect of condensed counterions is sim- 
ply to reduce uniformly the backbone charge density. In 
contrast to OSF result, the effects of condensation can- 
not be described in terms of a single renormalized pa- 
rameter q. The first term in eq 1111 arising from the net 
charge repulsion, tends to stiffen the chain, while the 
second term, originating from charge fluctuations, soft- 
ens the chain and is responsible for the chain collapse at 
low temperatures 0, 0, • 

In the limit of small nb, the persistence length is domi- 
nated by large-scale charge fluctuations, corresponding to 



small-fc contributions. This leads to the following asymp- 
totic expression for the persistence length 



^OSF 



£ 2 A 2 



(1 + 2£A K (nb)) 2 16 l + 2£AK (nb) 



(12) 



Note that the charge correlation contribution, i.e., the 
second term, does not exhibit simple scaling behavior as 
a function of x = A£. For small x, it varies quadratically 
with x, but it crosses over to one that scales as x for large 
x. 

The OSF result is valid for A£ < 1 (Cf. Fig.QJ. Thus 
when the charge fluctuations are important, OSF result is 
not accurate and the fluctuation correction should be in- 
corporated properly, as will also be detailed later. Our re- 
sults in eq^l should be compared with those for polyam- 
pholyte chains reported before (25|. In the former case, 
we sum up all multipole expansions [32j, while in the 
latter case the free energy is expanded in p owers of q 
and A. Inclusion of all multipole terms jlOl llll ] is cru- 
cial in the present case. This is because the counterion- 
mediated attraction becomes important when the chain 
is highly charged, i.e., for large values of / p. In this case, 
the multipole expansion diverges 0, On the other 
hand, q and A can be adjusted independently in the case 
of polyampholytes, and there can be finite ranges of pa- 
rameters where the perturbative description is valid. The 
convergence of the expansion is ensured for the screened 
cases if A is small enough, as is the case for the prob- 
lem discussed in Ref [25|. More recently, Golestanian et 
al. [2(j and Ariel and Andelman j2?J have also consid- 
ered a similar problem of the bending rigidity of a stiff 
polyelectrolyte. In contrast to the theory presented here, 
they only consider single-loop cases. At a low salt limit, 
our attractive bending rigidity scales as k" 1 log(l/«;6), 
which is in agreement with their result. On the other 
hand, our result in eq^] differs from the one in Ref. p^ . 
which is based on a SCL model of condensed counte- 
rions. If our result includes charge correlations driven 
by long-wavelength fluctuations, then the SCL approach 
captures short-wavelength charge correlations. Note that 
these two kinds of charge correlations should contribute 
to £ e iec- More precisely the competition between the two 
will essentially determine the leading behavior of the per- 
sistence length. In Appendix A, we provide a heuristic 
argument to study this competition. We find that the 
persistence length can be mainly determined by long- 
wavelength correlations as captured by our approach, 
even when the charge correlation between counterions 
is dominated by short-wavelength correlations. This is 
because small bending is more efficiently felt by long- 
wavelength charge correlations. In other words, the mode 
of charge correlations that dominates £ e i ec also depends 
on the mode of chain deformation, not just on counterion 
ordering. 

The electrostatic free energy in eq0]depends on the av- 
erage number of the condensed counterions per monomer, 
f a (a = 1,2). To solve for f a self-consistently, we equate 
the chemical potentials of condensed and free ions [20| . 
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In the following derivation, the subscript a — 1 and 2 
refer to monovalent and multivalent ions, respectively. 
The chemical potential of free ions is mainly associated 
with the configurational entropy of mixing and is approx- 
imately given by n^ ee /ksT ~ logn Q t>o, where n a is the 
bulk ion concentration and vq is the volume of ions. The 
chemical potential of condensed counterions arises from 
the attraction of counterions to PE chains and the en- 
tropic penalty for confinement |33|: 



..cond ..flue 

1 - -Z a ZqK {bK)+logf a 



k B T 



k B T 



, a =1,2, (13) 



where /i£ iuc is the charge fluctuation contribution to the 
chemical potential given by 



/■'-; 



flue 



k R T 



d \ f, f°° dkr ( _ 



xKoib^K 2 +k 2 )) - 2A^K (bVn 2 + k 2 )] 1.(14) 



Note that the term in {. . .} is a standard charge fluc- 
tuation free energy of an one-dimensional Debye-Huckel 
system immersed in a three dimensional ionic fluid |l0j| . 
Strictly speaking, our expression for [f° nd is valid for a 
rodlike chain conformation. In principle, the effect of 
chain deformation can be included; \i c ° nd can also be ex- 
panded in powers of 1/R. The leading-order behavior of 
£ e iec, however, is not influenced by the chain deformation 
correction to since it only leads to a subleading 

correction to £ e iec- The equilibrium value of f a can be 
fixed by requiring ^J ee — \i c ° nd . Note that can be 

important when Z is large and k _1 is finite |34|. 

To demonstrate the potency of multivalent counterions 
in softening PE chains, we have solved for f a and the elec- 
trostatic persistence length £ e i ec simultaneously. We plot 
leiec as a function of fo in Fig.^ We have chosen the pa- 
rameters n\ = ImM, b = 1.7A, r c = 2A, and £b = 7.lA 
(corresponding to T = 300K in water for which e = 80). 
First consider the case for which counterions are mono- 
valent (i.e., tt.2 = 0) as described by the thin curves. In 
this case, £ e i ec changes non-monotonically as fo increases. 
Our calculation should be compared with the correspond- 
ing OSF result £osf, which varies monotonically. Note 
that the OSF curve does not vary quadratically with fo 
beyond fo — 0.2. This is because condensed counterions 
start to renormalize the backbone charge density beyond 
fo ~ 0.2. The discrepancy between OSF result and ours 
may seem surprising. Our theory predicts that the chain 
size increases as the strength of Coulomb interaction in- 
creases, up to /o — 0.2. Beyond this, the chain shrinks its 
size with the increasing backbone charge fraction. This 
puzzling behavior can be understood in terms of the com- 
petition between the net charge repulsion and counterion- 
mediated attractions. When fo < 0.2, the long-ranged 
repulsion is dominant and the size of chain grows approxi- 
mately quadratically with fo, consistent with OSF result. 
The OSF behavior crosses over to the attractive regime 
where the charge fluctuations start to shrink the chain 



size. For sufficiently large fo, £ e iec can become negative. 
In other words, the bending rigidity of highly charged PE 
chains is smaller than the corresponding non- ionic chains. 
The reduction in the bending rigidity is attributed to the 
charge correlation effect that becomes dominant over the 
respulsive contribution due to a finite excess charge. 
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FIG. 1: The electrostatic persistence length £ e i ec (the solid 
lines) as a function of fo (the fraction of charged monomers). 
We have chosen the parameters b = 1.7A, £b = 7.lA, and 
ni = ImM. Due to the competition between the net charge 
repulsion and charge correlations, l e uc does not vary mono- 
tonically. In the absence of multivalent ions (712 = 0), £ e iec 
grows monotonically with fo up to fo = 0.2, consistent with 
OSF result (the thin dotted line). Beyond this, however, 
charge correlations become important and tend to shrink the 
persistence length. In the presence of a small concentration of 
O.lmM of trivalent counterions (Z = 3), charge correlation ef- 
fects are dominant and tend to reduce £ e iec beyond fo — 0.05. 
Over a wide range of fo (fo > 0.1), £ e iec is negative. Due to 
preferential adsorption of multivalent counterions, even OSF 
result (the thick dotted line) deviates significantly from the 
corresponding monovalent case (the thin-dotted line) 

The effect of counterion condensation is far more pro- 
nounced in the presence of a small concentration of 
O.lmM of trivalent counterions (Z = 3). Both OSF re- 
sult (the bold dotted line) and our result (the bold solid 
line) start to deviate from the monovalent case (corre- 
sponding to the thin curves) beyond fo ~ 0.05. In other 
words, the presence of O.lmM of trivalent counterions 
is more influential on the bending rigidity of PE chains 
than that of ImM of monovalent ions. The result for 
£ e iec in this case is strikingly different from the corre- 
sponding monovalent case; £ e iec is much smaller than in 
the corresponding monovalent case as long as fo > 0.05; 
it is also negative over a wider range of fo. Our results 
clearly demonstrate the efficiency of multivalent coun- 
terions in shorterning the persistence length of the PE 
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FIG. 2: The electrostatic persistence length (described by 
the solid line) as a function of ni (the concentration of Z- 
valent counterions) for b — 1.7A, Is = 7.1A, n\ = ImM, and 
Z = 3 (trivalent salts). The persistence length £ e i ec varies 
non-monotonically with ri2 for the all cases of fo = 1, 0.75, 0.5 
and takes on a minimum at a certain value of ri2 in the range 
O.OlmM < 712 < O.lmM; only a small concentration is needed 
to soften PEs. The efficiency of trivalent counterions in soft- 
ening the PE chain arises from the interplay between prefer- 
ential adsorption of multivalent counterions and charge fluc- 
tuations in determining £ e i ec . Also note that there is an op- 
timal value of 7i2 at which £ e i ec has the minimum value. At 
zero concentration of PE chains as in this figure, only a tiny 
concentration of multivalent counterions is needed to almost 
completely neutralize the backbone charges of the PE chain. 
Beyond the optimal value, multivalent ions start to contribute 
toward screening the charge correlation effect, leading to the 
non-monotonic variation of £ e iec with nz. 



chain, thus causing the collapse of the PE chain. The 
efficiency of Z = 3 in reducing £ e i ec can be attributed 
to the interplay between preferential adsorption of mul- 
tivalent counterions (over monovalent ones) and charge 
fluctuation effects in determining the persistence length; 
highly charged macrions can preferentially bind multi- 
valent counterions even when the bulk concentration is 
dominated by monovalent ions. Thus multivalent coun- 
terions are far more efficient on a molar basis in neutral- 
izing the macroion charge The OSF result in the 
presence of a tiny concentration of multivalent counteri- 
ons can strongly deviate from that for monovalent cases. 
Additionally, the strength of charge fluctuations grows 
approximately linearly with the valency of counterions, 
if preferential adsorption is assumed. When these two 
effects are combined, multivalent counterions can dra- 
matically change the bending rigidity of PE chains as 
demonstrated in Fig. ^ 

To further demonstrate the dramatic effects of multiva- 



lent counterions, we have estimated the electrostatic per- 
sistence length as a function of the concentration of mul- 
tivalent counterions n 2 . We have chosen the parameters 
Z = 3, m = ImM, b = 1.7A, r c = 2k, and £ B = 7.lA. 
We plot the electrostatic persistence length £ e i ec in Fig. [21 
as a function of n 2 for a few different values of fo- Our re- 
sults in the figure are intriguing. For all these cases (i.e., 
fo = 1, 0.75, 0.5), £ e iec varies non-monotonically with n 2 . 
It becomes more negative as n 2 increases from zero up 
to a certain value of n 2 . Beyond this, it decreases in 
magnitude with increasing n 2 - This implies that there is 
an optimal concentration at which l e i ec is most negative. 
Somehow this optimal value is not sensitive to /o; it falls 
in the range O.OlmM < n 2 < O.lmM and hence only a 
small concentration is needed to soften PE chains. This 
clearly demonstrates the efficiency of multivalent coun- 
terions in softening oppositely charged PE chains. For 
the entire range of n 2 adopted in Fig. [21 £ e iec does not 
show a simple scaling law. When n 2 is somewhat smaller 
than the optimal value, £ e iec has weak dependence on n 2 
as shown in the figure. We find that, for all cases (i.e., 
fo = 1,0.75,0.5), £ e i ec assumes the simple scaling form: 
^eiec ~ c + d log n 2 , where c and d are ^-independent 
negative numbers and the concentration n 2 is in mM. 
As evidenced in Fig. [21 d (i-e., the slope of the curves) 
is insensitive to fo] it has been estimated to be close to 
— 0.2A for all three cases. In contrast, c changes linearly 
with fo and is more negative for larger fo; for the chosen 
parameters, c has shown to be approximately given by 
c w — 10/o (A). It should be emphasized that both c and 
d are non-universal constants that can depend on k and 
Z for fixed £b and b. Also note that this simple scaling 
behavior can be invalidated by the competitive binding 
on the one hand and the presence of the opposing effects 
in determining £ e i ec (cf. eg 1111) on the other. 

The efficiency of multivalent counterions as evidenced 
by our results in Fig. [3 can also be understood in terms 
of competitive binding and charge fluctuations. When 
n 2 = 0, only the monovalent counterions bind to the PE 
chain. When n 2 ^ 0, however, the PE chain can preferen- 
tially bind multivalent counterions, replacing monovalent 
counterions in its vicinity, thus softening the PE chain, as 
can be seen by solving the linear equations in ea ll3l simul- 
taneously [3(|. For sufficiently highly charged PE chains 
(fo ~ 1) at zero concentration, only a tiny concentration 
of multivalent counterions is needed to replace monova- 
lent counterions. Note that the optimal value of n 2 at 
which £ e i ec has the minimum value depends not only on 
n\ but also on the concentration of PE chains. Beyond 
the optimal value, multivalent counterions start to con- 
tribute toward screening, weakening the correlation effect 
and thus leading to the non-monotonic variation of £ e iec 
with n 2 . Our result in Fig.[21corresponds to PE chains at 
zero concentration and thus a direct comparison of our 
result to experimantal data for PE chains at non-zero 
concentration should be made with due caution. Finally, 
computer simulations have explicitly provided evidence 
for the efficiency of multivalent counterions in inducing 
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collapse of strongly charged polyelectrolytes 01 • 



C. Multi- loops 

So far we have restricted ourselves to the single-loop 
case. Due to the breakdown of the pairwise additivity of 
charge-fluctuation interactions 1 1 OL 1 1 1| . it is important to 
discuss the effects of inter-loop coupling on the bending 
rigidity. Based on our common intuition, we expect the 
inter-loop interactions to stiffen the condensate. This is 
because the net charge repulsion is longer ranged than 
the attraction This conjecture is, however, based on 
the pairwise additivity of the electrostatic interactions 
between different loops and can easily be invalidated un- 
less the charge fluctuations are small enough fiol Ull l3l| . 
For the case of parallel rods, it has been shown that the 
pairwise additivity is valid only when the interaction is 
repulsive [T(J, UM- By the same token, it is crucial to 
include non-pairwise additive interactions or many-body 
effects 37] in the description of the bending rigidity of 
PE condensates. Because the solution of eq^Jin general 
cases is difficult, we invoke physically motivated simplifi- 
cation. To understand the physical consequences of non- 
pairwise additive interactions on the bending rigidity, we 
make some approximations. First we do not explicitely 
include excluded volume repulsion between loops. To 
keep different loops from approaching arbitrarily close to 
each other, we assume that loops in the condensate are 
arranged on a square lattice with lattice constant a. The 
lattice constant a can be considered as an equilibrium 
distance between two neighboring loops, which is deter- 
mined by the balance of a few competing effects: charge 
correlation attractions, excluded volume repulsions, etc. 
We also assume that all loops in the bundle remain par- 
allel with each other, even when they bend, and that Aj 
and cjj are independent of j and taken to be equal to A 
and q, respectively. With this simplification, we can take 
advantage of this periodicity to recast the problem using 
discrete Fourier transform. The resulting bending free 
energy is given by 



q 



s/N-l 
■ - M f ~~ E 

2£A [ K Q )\ hh=o 
dk 2A£A' (kL,/c) 



27ri + 2A^ (k_L,fe) 



J dsDjJs) 



cos(k_L ■ j±) 



x cosfcs, 



(15) 



where kj_ is the wave vector conjugate to ajj^, where 
j± = (ji, h) and ji, j 2 = 0.1, . . . , \/N - 1, and is given 
by k ± = -^L(m,n 2 ) with n u n 2 = 0, 1, . . . , VN - 1. 
We denote the discrete Fourier transform of a function 
9(a\j±\VK 2 +k 2 ) as 



<?0,kj 



£ 



Finally, the function °g is simply g(k = 0,k± =0). 

First note that the bending free energy of iV-loop con- 
densates is not simply N times that of a single loop. 
Our theory suggests that each loop in a bundle is fur- 
ther softened due to the inter-loop coupling, as implied 
by eq ED We expect that the repulsive bending rigid- 
ity per loop decreases with N. This is because coupling 
between charges on different loops enhances screening in 
electrostatic repulsions. Similarly, the attractive interac- 
tion becomes stronger as iV increases, due to enhanced 
charge correlations as will be detailed later. More pre- 
cisely, inter-loop correlations enhance intra-loop correla- 
tions. A similar issue for the case of rigid PEs is discussed 
by Ha and Liu |3J| . The implication of eq ED is most 
striking in the limit of nb — > 0: For sufficiently large N, 
the repulsive term in eq 1151 scales as 1 /N and becomes 
vanishingly small as N — > oo. Thus, only the attraction 
can modify the bending rigidity in this case. Recently, 
Ha and Liu |3ll| have considered two interacting bundles 
of randomly-charged rods and shown that, as Kb — ► 0, 
the interfacial interaction comes from charge-fluctuation 
attractions only. This is analogous to the vanishingly 
small contribution of the repulsion to the bending rigid- 
ity in the present case. For the screened case of n ^ 0, 
on the other hand, the repulsion also contributes to the 
bending rigidity. 

For sufficiently large N, we can show that the con- 
densate has a well-defined bulk bending free energy: 
A^eiec ~ Vibuik where V is the condensate volume and 
(■bulk is iV-independent persistence length. To establish 
this, we examine the behavior of Kq, Ko(k, kj_), and 
D(k,k±) for large N. For large a;, Kq(x) ~ l/y/xe~ x and 
thus °Kq, °D, Ko(k, kj_) approach constants as N — > oo. 
Similarly, it can be shown that °D, D(k, kjj also become 
constants as N — > oo. This is trivially true for finite k. 
As k — > 0, it suffices to establish this asymptotic behav- 
ior for the attractive term only, since the repulsive term 
becomes dominated in this limit. As n — > 0, Dij(s) in 
eq [3 vanishes at least as fast as 1/-R^ for large \i — j\. 
Thus the descrete Fourier transform is convergent. Fur- 
thermore, for large N, we can replace the sum J^k^ can 
be replaced by N (f^jr, up to a correction of or- 

der unity. From these arguments, we conclude that the 
bending free enegy per monomer increases linearly with 
N for sufficiently large N. This enables us to write the to- 
tal bending free energy as AF e i ec /k B T = \l bu ikV '/ R 2 a 2 , 
where (bulk is A^-independent persistence length per loop. 
An asymptotic behavor oiibuik can be obtained replacing 
the summations by integrals; we found that, for kD <C 1, 
hulk varies as 



(-bulk 



Iog(«-VD)' 



(17) 



g(a\h\ V k 2 + fc 2 )cosk_L ■ (16) 



where D is the crosssectional diameter of the bundle, 
D = \fNa. This should be compared with the corre- 
sponding result for a single loop case; the PE chains are 
further softened by inter-loop couplings by the factor of 
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log(6/D). Park et al. |23 have used a similar expression 
for the bending free energy (Cf. their eq 1) without the 
benefit of derivation. Here, we have shown that "many- 
body effects", i.e., loop-loop couplings, lead to the bend- 
ing free energy of PE condensates which grows linearly 
with N for N ^> 1, providing a quantitative basis for the 
previous work of Park et al. [23|. 

III. CONCLUSIONS 

In conclusion, we have studied the effect of charge 
fluctuations on the stiffness of highly charged stiff PE 
chains. In particular, we have focused on the interplay 
between competitive binding of counterions to PE chains 
and charge correlations between condensed counterions 
in softening the PE chains. We have shown that the 
bending rigidity of the PE chain, in the presence of coun- 
terion condensation, cannot be captured by a simple scal- 
ing behavior. The existence of multiple distinct regimes 
characterizes the conformations of highly charged PEs. 
This is a consequence of the simultaneous presence of 
competing interactions, which tends to invalidate a sim- 
ple scaling analysis. Our theory also illustrates the signif- 
icance of non-pairwise additivity of counterion-mediated 
interactions in softening PE chains; we have shown that 
the inter-loop coupling enhances softening, resulting in a 
well-defined bulk bending free energy. 

Acknowledgements: This work was supported in part 
by the National Science Foundations through grant num- 
ber CHE0209340 (DT) and by the Natural Science and 
Engineering Research Council of Canada (BYH). We are 
grateful to D. Andelmann for helpful comments. 



Appendix A 

In this appendix we argure that the Debye-Huckcl 
(DH) approach may be valid over a much wider parame- 
ter space than implied by a simple thermodynamic con- 
sideration in Ref. [23. In Ref. [2||, condensed counteri- 
ons are considered as forming a strongly correlated liquid 
(SCL) confined to the surface of their binding PE chains. 
It was then argued that the negative, charge correlation 
contribution to the persistence length is dominated by 
these strong charged correlations of the SCL. In this 
appendix, we argue the thermodynamic argument un- 
derestimates the importance of long-wavelength charge 
fluctuations as compared to SCL correlations, since it 
ignores the coupling of charge correlations to bending. 
Recently, it has been shown that, at high temperatures, 
long-wavelength (LW) charge fluctuations dominate the 
free energy, while short- wavele ngt h (SW) fluctuations are 
dominant at low temperatures |39j . Even below the freez- 
ing temperature, there exists a LW contribution to the 
free energy. Interestingly the major contribution to the 
electrostatic bending rigidity of a polyelectrolyte chain 
can arise from LW charge fluctuations even when the 



free energy is dominated by SW fluctuations, as long as 
the chain is near the rod limit. This is because small 
bending is more effectively felt by the LW charge fluctu- 
ations. In this case the DH approach ought to be a good 
approximation. 

To focus on the essential physics of this issue, we here 
consider a single-loop case where all backbone charges are 
neutralized by counterions. Near the rod limit, we can 
write the interaction Hamiltonian as follows: TL — 7i r od + 
ATI where 7i r od corresponds to the rodlike conformation 
and AH. is the change in the interaction due to bending. 
The free energy cost due to bending is 

AT = {AH) = \Ys D ^ s ') (<?( s M s ')> + ° ( R ~ 4 ) > 

ss' 

where (• • • ) is an average with respect to the Boltzmann 
factor e - H raA/k B T and £,( Sj s ') is 

The charge correlation contribution to the persistence 
length can be read off from this expression: 

£oorr = ^^0(3,3') (q(s)q(s')). (A3) 

ss' 

Now the computation of the persistence length reduces 
to the computation of the charge correlation function 
(q(s)q(s')) . So far our formalism is exact up to 1/R 2 . 

If we use a Debye-Huckcl (DH) approximation, then 
(q(s)q(s')) is simply 

(q(s)q(s')) DH = Q- 1 (s,s% (AA) 

where Q is a matrix whose matrix element Q(s, s') is 

Q(s, s') = A~ 1 5 ss - + £ B ° ^ '— . (A5) 

\s - s'\ 

As a result, the persistence length in eq A3 reduces to the 
one in ea !12l In the limit of n — > 0, the DH approximation 
leads to 

^ (CA) 2 

" ~ 16 l + 2£Alog(l/«6)- {M) 

The persistence length in this expression grows linearly 
with the Debye screening length. 

At low temperatures, the backbone charges together 
with the condensed counterions tend to be strongly 
correlated. To simplify the problem we assume that 
/o = 1 and that counterions are divalent (Z = 2) 
and form an ionic crystal such that the charge distri- 
bution along the chain is represented by a ground state 

(H 1 1 1 — ). The resulting charge correlation can 

then be approximated by the following oscillatory func- 
tion: (q(s)q(s')) osc ~ cos(7r|s — s'\/b). If we use this 
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expression in eq A3, then we have 



get a simple criterion for the DH approach to be valid: 



6^ 

s— 1 



In the limit of Kb — > 0, 



(s 2 + ks 3 ) cos ( — 



approaches 

£ = - l JL 

26' 



(A7) 



(AS) 



This result can be compared to the one based on the SCL 
model of polyelectolytes HJ, Iscl = -Z 2 £ B /96. For 
Z = 2, this leads to ^ S cl = ~W 24 - Note that £ S cl 
is only slightly different from the one in eq ^48. This 
difference can be attributed to the fact that we treat the 
charge fluctuations differently from the SCL approach. 
Note that Idh is much larger than £scl for small nb. 
This is becase the LW charge fluctuations are much more 
effectively felt by small bending in this case. 

However, neither of £dh nor £osc solely describes the 
bending rigidity of polyelectrolyte chains for a wide range 
of T or f . This is because lose suppresses fluctuations 
while £uh does not accurately capture strong charge fluc- 
tuations. In general both LW and SW charge flcutuations 
contribute to the persistence length. To study the com- 
petition between the two, we use a linear response theory 
and assume that the free energy of the charge fluctuations 
on a rod can be written in the Fourier space as 



^rod 



££>(*)S -1 (*)*9(- 



■k), 



(A9) 



where S(k) is the charge structure factor. The probabil- 
ity of the fluctuation Sq(k) is proportional to the factor 
exp [— S~ 1 (k)] . At high temperatures, S(k) = SDH(k) — 
Q^ 1 (k) while at low temperatures S(k) = Sosc(k) — 
<5fc iL . At intermediate temperatures, we assume that 
S(k) ~ A(T)S DH {k) + B(T)S sc{k). Within this ap- 
proximation, the charge correlation persistence length is 
given by 



- A(T)£ DH - B(T)l OS c 



(AW) 



Suppose the free energy is dominated by short- 
wavelength fluctuations, i.e., A(T) -C B(T). Even in this 
case, however, the persistence length can be dominated 
by long- wavelength fluctuations, i.e., £ CO rr ~ A(T)£dh 
in the limit of Kb — > 0. This is because small bending is 
much more effectively felt by the long-wavelength fluctu- 
ations in this limit. Unfortunately, the precise form of 
A(T) and B(T) for a wide range of T are not known to 
date. Nevertheless it is clear that A(T) w 1 is a good 
approximation at high temperatures, while B(T) — ► 1 as 
T — > 0. Near and at the crossover region between £dh 
and £osCi however, A(T) and B(T) may deviate from 
this asymptotic value, i.e., 1. If we use this, we would 



log(l/«6) 



> 1. 



(All) 



Note that this can easily be satisfied at a low salt limit. 
This implies that the electrostatic persistence length can 
be mainly determined by Idh, even when A(T) w B(T). 
As a result, our DH approach is valid for a much wider pa- 
rameter space than implied by Ref. ,29j which suppressed 
the interplay between chain deformation and charge cor- 
relations in determining the mode of charge correlations 
that dominates the persistence length. 



IV. APPENDIX B 

In this appendix, we outline the derivation of our major 
result in eq 1101 In this appendix we use D to denote 
Q Qo- First consider the first term of ea 117)1 To simplify 
this term, consider 



/^[Qo^DQo 1 ]..,^ 

ss' 



(SI) 



In this equation and in what follows, (A ss i ) ■-, = 
A- ss >jf, where A is a block matrix such as Q. It 
then follows that [(AB) ss ,] jf = E^^B^^ = 
Yljii Y] 3 n Ajfi as itBftji s ti a t, where B is also a block ma- 
trix of the same rank. Note that (Qo ) JJ " SS " IS a func- 
tion of \s — s"\. It is easy to show that 



°M7i„ (B2) 



where °A4 is the same matrix defined in eq 8. Here, the 
second step follows from 



E^oV 



ii" 



l^l \ A-i 4 



A- 1 + V 
^ A (1 - Ay + A 2 V 



ii" 

2 



-l ii" 



- [A (1 - A2£K + A 2 (2^ ) 2 + ...)]. fl 



where V jfss < = Zb '^^.,'^ and (K ) J 
Ko(KRjjt). If we use this in eq Bl, then we have 



(B3) 



X- 1 ]T Dm^m- 1 



(S4) 



J ii' 
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This , if combined with eq 6, essentially leads to the first 
term of eq 10 in the continuum limit. 

To simplify the second term of eq 10, consider 



J (Qq jjss 



(B5) 



jj 



Here (Qo 1 ) s 

, can be expanded as follows: 
1 \ 



(Qo 1 )^ 



Af + V 

\ J / ss' 

= AjUss' - AjV ss > + A 2 j^2 Vss"V s » s < + ... 

(B6) 

Using this in eq B5, we have 

J=Aj D ss -J2[ A jV ss >D ss , 

s' ^ 

— V SS " V S " s ' D s t s 

s' s" 

V ss" V 

s's" s <" '^ii 

(B7) 

Note that the first-order term in Aj in this equation will 
cancel the term containing A in the second term of eq 



10. As a result, the second term of eq 10 becomes 
2nd term = -\^A ( A 3 V ™' 

— Aj V ss" V s" s' D s > s 

s" 

+ Aj V ss"V S"S"'V S"'S'D 

s" s'" 33 

(B8) 

In the Fourier space k, conjugate to s — s' , this equa- 
tion can be diagonalizcd with respect to s and s' and 
resummed: 

2nd term = -yp,|| ^ (aj V(k)D(k) 

3 

- MV(k)V(k)D(k) 



+ A]V{k)V{k)V{k)D{k) 



M V- x 



rffc A j V(fc)£>(fc) 



2tt l + AjV{k) 



33 

where f(k) is a Fourier transform of f(s): f{k) — 
12 s /( s ) cos This essentially leads to the second term 
of eq 10. 
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In other words, (. . .) reduces to (. . .) if the electrostatic 
interactions among counterions are turned off. In this 
case, we can consider condensed counterions as forming 
an ideal gas in chemical contact with free counterions. As 
a result, charges at different sites are uncorrelated. If M( 
and M| are the total number on loop j of monovalent 
and multivalent counterions, respectively, then it follows 
from a simple thermodynamic consideration that 



(Ml 



Ah 



Ml 



(19) 



where a — 1,2. This gives rise to charge fluctuations 
along the chain; the total-charge fluctuation of condensed 
counterions is 



MA; 



E*2 



(Ml 



Mf- 



= ZlMl 



where Z\ = 1 and Z2 — Z. This leads to the variance in 
the charge per site charge used in Sec. II: Aj = fi + Z 2 fa- 
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